Screening of Reference miRNA of Different Early- and Late-Flowering Tree Peony Varieties

miRNA plays an important role in plant growth and development and in response to various stresses. Quantitative real-time PCR (qRT-PCR) technology is often used to detect the expression level of miRNAs and genes by comparing with reference genes. In order to screen out the optimal reference miRNAs in different tree peony varieties, the petals of 42 different early- and late-flowering tree peony varieties were used as experimental materials, and geNorm, NormFinder, Bestkeeper, and RefFinder software were used to evaluate the stability of 16 candidate reference miRNAs. The results showed that the average Ct values of all candidate reference miRNAs were between 15.34 ± 0.29 and 32.64 ± 0.38. The optimal number of reference miRNAs was four, which were PsPC-5p-19095, PsPC-3p-51259, PsmiR159a, and PsPC-3p-6660 in geNorm. The stability of PsPC-3p-6660 was the highest in the analysis results of NormFinder software. Among the analysis results of Bestkeeper software, PsMIR319-p5 has the highest stability. Among the results of comprehensive evaluation and analysis of several software using RefFinder, the candidate reference miRNA with the highest stability was PsPC-3p-6660. When PsPC-3p-6660 was used as the reference miRNA, the expression of PomiR171 and PomiR414 in response to different flowering times of tree peony was relatively stable in 42 tree peony varieties, indicating that PsPC-3p-6660 was stable and reliable. The results of this study provide a reference miRNA for studying the expression changes of miRNA in different tree peony varieties and further exploring the regulatory mechanism of miRNA in different peony varieties.


Introduction
Tree peony (Paeonia suffruticosa Andrews.) is a perennial woody plant, which belongs to Paeoniaceae, Paeonia, section Moutan. As a kind of multi-purpose plant with economic value, tree peony has high ornamental value because of its rich varieties, gorgeous colors, and diverse flowers [1,2]. In addition, it also has high nutritional and medicinal values because its root can be used as medicine and the seeds can be used to extract oil, which is beneficial to human health [3][4][5][6]. The natural flowering period of this tree is around April, which is a one-flower-a-year plant, and a few varieties have the phenomenon of second flowering [7,8]. Among them, middle-flowering tree peony varieties dominate, and early-flowering and late-flowering varieties are few. The characteristics of a relatively short and concentrated flowering period seriously affect the ornamental value of tree peony and limit the development of tree peony industry. Therefore, understanding the internal molecular mechanism of tree peony flowering may be helpful to discover new genes related to flowering traits, which can be used to regulate the flowering time of tree peony, prolong the group flowering period, and increase economic benefits.
At present, the research on tree peony mainly focuses on the construction of tissue culture rapid propagation system [9], cultivation and domestication [10], genetic diversity [11],

miRNAs Quality Analysis
The miRNAs of 42 different early-and late-flowering tree peony varieties were extracted from the petals. The integrity of miRNAs was detected by 3% agarose gel electrophoresis, and the concentration and purity of miRNAs were detected by NanoDrop 1000 spectrophotometer (Implen, Munchen, Germany). The miRNA of each sample had a single bright and clear band at about 100 bp, indicating that the integrity of miRNAs was good. The NanoDrop 1000 spectrophotometer detection results showed that the A 260/280 nm of each miRNA sample was 1.8-2.1, indicating that the purity of miRNAs was good, which met the requirements of subsequent tests.

Primers Specificity Analysis of 16 Candidate miRNAs from Different Tree Peony Varieties
The primer melting curves of 16 candidate reference miRNAs were analyzed by qRT-PCR, and the results showed that the melting curves of all primers showed a single peak, indicating that the primers of the 16 candidate reference miRNAs had good specificity and could be used for the evaluation of the reference miRNAs ( Figure 1).

miRNAs Quality Analysis
The miRNAs of 42 different early-and late-flowering tree peony varieties were extracted from the petals. The integrity of miRNAs was detected by 3% agarose gel electrophoresis, and the concentration and purity of miRNAs were detected by NanoDrop 1000 spectrophotometer (Implen, Munchen, Germany). The miRNA of each sample had a single bright and clear band at about 100 bp, indicating that the integrity of miRNAs was good. The NanoDrop 1000 spectrophotometer detection results showed that the A260/280 nm of each miRNA sample was 1.8-2.1, indicating that the purity of miRNAs was good, which met the requirements of subsequent tests.

Primers Specificity Analysis of 16 Candidate miRNAs from Different Tree Peony Varieties
The primer melting curves of 16 candidate reference miRNAs were analyzed by qRT-PCR, and the results showed that the melting curves of all primers showed a single peak, indicating that the primers of the 16 candidate reference miRNAs had good specificity and could be used for the evaluation of the reference miRNAs ( Figure 1).

Figure 1.
Primer melting curves analysis of 16 candidate reference miRNAs in different tree peony varieties. A single peak indicates that the primer specificity is good.

Ct Value Analysis of 16 Candidate miRNAs of Different Tree Peony Varieties
The expression levels of miRNAs in 16   Primer melting curves analysis of 16 candidate reference miRNAs in different tree peony varieties. A single peak indicates that the primer specificity is good.

The Expression Stability of 16 Candidate Reference miRNAs Analyzed by geNorm
The geNorm software mainly analyzes the expression stability of candidate reference genes through the expression stability M. The larger the M value, the more unstable the gene, and the smaller the M value, the more stable the gene. The analysis results showed that the stability order from high to low of 16 candidate reference miRNAs was PsPC-5p- 5p-19095 and PsPC-3p-51259 had the best stability, and the most unstable expression was PsmiR171k-3p (M = 1.09) ( Figure 3). Pairwise variations can be used to determine the optimal number of reference genes. The results of geNorm histogram analysis showed that V4/5 = 0.13 < 0.15, indicating that the optimal number of references for miRNA quantitative expression was 4, which were PsPC-5p-19095, PsPC-3p-51259, PsmiR159a, and PsPC-3p-6660 ( Figure 4).

The Expression Stability of 16 Candidate Reference miRNAs Analyzed by geNorm
The geNorm software mainly analyzes the expression stability of candidate reference genes through the expression stability M. The larger the M value, the more unstable the gene, and the smaller the M value, the more stable the gene. The analysis results showed that the stability order from high to low of 16 candidate reference miRNAs was PsPC-5p- 5p-19095 and PsPC-3p-51259 had the best stability, and the most unstable expression was PsmiR171k-3p (M = 1.09) ( Figure 3). Pairwise variations can be used to determine the optimal number of reference genes. The results of geNorm histogram analysis showed that V4/5 = 0.13 < 0.15, indicating that the optimal number of references for miRNA quantitative expression was 4, which were PsPC-5p-19095, PsPC-3p-51259, PsmiR159a, and PsPC-3p-6660 ( Figure 4). ' in the box is the outliers, the horizontal line in the box is the median, and the upper and lower side lines of the box are the upper four digits and the lower four digits, respectively.

The Expression Stability of 16 Candidate Reference miRNAs Analyzed by geNorm
The geNorm software mainly analyzes the expression stability of candidate reference genes through the expression stability M. The larger the M value, the more unstable the gene, and the smaller the M value, the more stable the gene. The analysis results showed that the stability order from high to low of 16 candidate reference miRNAs was PsPC-5p indicating that PsPC-5p-19095 and PsPC-3p-51259 had the best stability, and the most unstable expression was PsmiR171k-3p (M = 1.09) ( Figure 3). Pairwise variations can be used to determine the optimal number of reference genes. The results of geNorm histogram analysis showed that V 4/5 = 0.13 < 0.15, indicating that the optimal number of references for miRNA quantitative expression was 4, which were PsPC-5p-19095, PsPC-3p-51259, PsmiR159a, and PsPC-3p-6660 ( Figure 4).

The Expression Stability of 16 Candidate Reference miRNAs Analyzed by NormFinder
NormFinder can rank the stability of candidate reference genes according to the stability value S of gene expression. NormFinder can not only compare differences between candidate genes but also calculate variations between sample groups. The smaller the S value, the more stable the gene, otherwise the worse the stability of the reference gene. According to the evaluation results of NormFinder software, the most stable expression of the 16 candidate reference miRNAs was PsPC-3p-6660 (S = 0.37), followed by PsPC-5p-19095 (S = 0.39), and the most unstable was PsmiR171k-3p (S = 0.79) ( Table 1). The expression stability of 16 candidate reference genes in different early-and late-flowering tree peony varieties was ranked from high to low as follows: PsPC-3p-6660, PsPC-5p-19095, PsMIR319-p5, PsPC-3p-51259, PsPC-5p-9292, PsMIR11609-p5, PsPC-3p-18408, PsmiR159a, PsPC-3p-23386, PsPC-3p-70893, PsPC-3p-13662, U6, PsmiR858-3p, PsmiR11607, PsPC-3p-15676, and PsmiR171k-3p. Because NormFinder software can only screen out one best reference gene, PsPC-3p-6660 was selected as the best reference base in the analysis using NormFinder software. The results were similar to geNorm software but with slight differences. The Bestkeeper software mainly evaluates the stability of the reference genes by calculating the standard deviation (SD) and the covariance (CV). The smaller the SD value, the better the gene stability. The analysis results showed that the SD values of most candidate miRNAs were less than 1, among which PsMIR319-p5 had the best stability, followed by PsPC-33-6660. Only one candidate miRNA had an SD value greater than 1, namely, PsmiR171k-3p, indicating that it had the lowest stability (Table 2). The smaller the SD and CV values, the more stable the gene is.

The Expression Stability of 16 Candidate Reference miRNAs Analyzed by RefFinder
Since the analysis results of different software have certain differences, the RefFinder software can be used for the final comprehensive analysis. RefFinder analysis is a comprehensive evaluation and ranking of the stability of reference genes obtained by geNorm, NormFinder, and Bestkeeper software. The analysis results showed that there were great similarities between the results obtained by RefFinder and the other three software, but there were also some differences. Among the petals of 42 different early-and late-flowering tree peony varieties, PsPC-3p-6660 was the most stable and consistent with the results of NormFinder, followed by PsPC-5p-19095. PsmiR171k-3p was the most unstable miRNA, and the stability order from high to low is as follows: PsPC-3p-6660, PsPC-5p-19095, PsMIR319-p5, PsPC-3p-51259, PsPC-5p-9292, PsMIR11609-p5, PsmiR159a, PsPC-3p-23386, PsPC-3p-18408, PsPC-3p-70893, PsmiR11607, PsPC-3p-13662, PsmiR858-3p, U6, PsPC-3p-15676, and PsmiR171k-3p. There was consistency with the most unstable miRNA obtained by the other three software, indicating that the results were reliable (Table 3). PsPC-3p-18408 8.10 10 PsPC-3p-70893 8.91 PsmiR171k-3p 16.00 The smaller the Geomean of Ranking Values, the more stable the gene is.

Validation of Reference miRNAs
In order to verify the reliability of the selected reference miRNAs, we selected the most stable reference, PsPC-3p-6660, and the most unstable reference, PsmiR171k-3p, as reference miRNAs to analyze the expression patterns of PomiR171 and PomiR414 in response to different flowering times in different tree peony varieties. The results showed that when PsPC-3p-6660 was used as the reference miRNA, the relative expression levels of PomiR171 and PomiR414 in different tree peony varieties were relatively stable, and the average Ct values of 42 petal samples were 0.24 ± 0.03~2.91 ± 0.55 and 0.08 ± 0.03~2.06 ± 0.68, respectively. When PsmiR171k-3p was used as the reference miRNA, the relative expression levels of PomiR171 and PomiR414 in different tree peony varieties were quite different. The average Ct values of 42 petal samples were 0.60 ± 0.03~26.30 ± 2.49 and 0.56 ± 0.31~32.18 ± 3.40, respectively ( Figure 5). Since PomiR171 and PomiR414 are differentially expressed miR-NAs in response to different flowering stages of tree peony, when PsPC-3p-6660 with high stability was used as the reference miRNA, the results of relatively stable expression of PomiR171 and PomiR414 among the samples were more reliable. When PsmiR171k-3p with the lowest stability was used as the reference miRNA, the results of relatively stable expression of PomiR171 and PomiR414 among the samples were not reliable. Therefore, PsPC-3p-6660 is the most suitable reference miRNA of different early-and late-flowering tree peony varieties. PsmiR171k-3p with the lowest stability was used as the reference miRNA, the results of relatively stable expression of PomiR171 and PomiR414 among the samples were not reliable. Therefore, PsPC-3p-6660 is the most suitable reference miRNA of different early-and late-flowering tree peony varieties.

Discussion
Because the differences between species, materials, and tissues will affect the stability and reliability of quantitative results, it is necessary to select appropriate reference genes for data standardization during qRT-PCR. In this study, 42 early-and late-flowering tree peony varieties were used as materials, and the best reference miRNA was successfully screened out, which provided suitable reference miRNA for miRNA research of different tree peony varieties.

Discussion
Because the differences between species, materials, and tissues will affect the stability and reliability of quantitative results, it is necessary to select appropriate reference genes for data standardization during qRT-PCR. In this study, 42 early-and late-flowering tree peony varieties were used as materials, and the best reference miRNA was successfully screened out, which provided suitable reference miRNA for miRNA research of different tree peony varieties.
In general, the same reference does not have stability all the time. Under drought stress, the most suitable reference miRNAs in roots and leaves of Glycine max were miR156a and miR167a [35]. In Allium sativum, the most stable reference miRNA for different explants was AsmiR168a-5p and for different genotypes, the most stable reference miRNA was AsmiR159a-1 [36]. In Juglans regia, the most suitable reference miRNAs for flower buds at different differentiation stages were jre-miR394a, jre-miR159a, and jre-miR159c, and the most suitable reference miRNAs for leaf buds at different differentiation stages were 5.8S rRNA and jre-miRn3 [37]. The most stable reference miRNA combinations during seed development in Brassica napus were miR167-1_2, miR11-1, and miR159-1 [38].
In this study, the stability analysis of 16 candidate reference miRNAs in 42 different tree peony varieties showed that PsPC-5p-19095 and PsPC-3p-51259 had the highest stability in the geNorm software. PsPC-3p-6660 had the highest stability in the NormFinder software. PsMIR319-p5 has the highest stability in the Bestkeeper software. In the results of comprehensive evaluation and analysis using RefFinder, the candidate reference miRNA with the highest stability was PsPC-3p-6660, which was consistent with the results of NormFinder and slightly different from the analysis results of geNorm and Bestkeeper software. The difference between the analysis results of each software has also appeared in the previous research results [39]. The reason for the difference in the stability of each candidate reference miRNA may be due to the large number of candidate reference miRNAs selected. Secondly, the difference in mathematical algorithms between different software will also affect the stability ranking of the test results to a certain extent. The results of several software showed that the stability of PsmiR171k-3p was the lowest, indicating that the candidate miRNA was not suitable as the reference gene in different tree peony varieties.
At present, U6 is one of the most common reference genes, which is used as the reference gene in miRNA quantitative expression analysis of various plants, such as Vitis vinifera [40], Brassica oleracea [41], and Jatropha curcas [42]. Studies have shown that U6 is not stable in all cases [43]. U6 is a suitable reference miRNA for different tissues and stem tissues under drought stress in Hylocereus polyrhizus [44]. However, in the evaluation of the expression stability of walnut flower bud and leaf bud differentiation, tissue parts, and varieties, U6 had the worst stability and was not suitable as the reference miRNA. This phenomenon also exists in plants, such as wheat and longan [45]. In this study, U6 was used as a candidate reference miRNA, and its stability was low in the four software, which was not suitable as the reference miRNA for different tree peony varieties.  (Table 4).
qRT-PCR primers were designed using poly (A) by Primer Premier 5.0 software. Primers should avoid special structures such as dimers and hairpin structures to prevent adverse effects on test results. Then, they were sent to Sangon Biotech (Shanghai) Co., Ltd. (Shanghai, China) for synthesis (Table 5). Table 5. Primers for qRT-PCR.

Primers Specificity Analysis
Using cDNA as a template, qRT-PCR primers were used for PCR amplification using a 2× PCR Taq

Expression Stability Analysis of 16 Candidate Reference miRNAs
The stability of 16 candidate reference miRNAs was evaluated and ranked using geNorm, NormFinder, and Bestkeeper software, respectively. And then the online tool RefFinder (http://blooge.cn/RefFinder/, accessed on 13 February 2023) [46,47] was used to comprehensively evaluate and analyze the results obtained by the above three software. Finally, the most suitable reference miRNA among different tree peony varieties was screened.

Validation of Candidate Reference miRNAs
The most stable and the unstable miRNAs were selected as the reference miRNAs, and the expression patterns of PomiR171 and PomiR414 (CNGBdb, accession number CNP0002984) in response to different flowering times of tree peony in different tree peony varieties were analyzed. Three technique replicates were set for each sample, and the qRT-PCR primers were given in Table 6. The expression levels of the miRNA were calculated using the formula of 2 −∆∆Ct , and statistical analysis and mapping were performed using SPSS 22.1 (one-way ANOVA, Duncan's test, p < 0.05) and Origin 2018.

Conclusions
Using petals of 42 different early-and late-flowering tree peony varieties as experimental materials, the expression stability of 16 candidate reference genes was evaluated and analyzed by geNorm, NormFinder, Bestkeeper, and RefFinder software using qRT-PCR. The results showed that the average Ct values of all candidate reference miRNAs were between 15.34 ± 0.29 and 32.64 ± 0.38. In the geNorm software, PsPC-5p-19095 and PsPC-3p-51259 had the highest stability, and the optimal number of reference miR-NAs was four, which were PsPC-5p-19095, PsPC-3p-51259, PsmiR159a, and PsPC-3p-6660. The stability of PsPC-3p-6660 was the highest in the analysis results of the NormFinder software. Among the analysis results of the Bestkeeper software, PsMIR319-p5 had the highest stability. Among the results of comprehensive evaluation and analysis of several software using RefFinder, the candidate reference miRNA with the highest stability was PsPC-3p-6660. Therefore, PsPC-3p-6660 can be used as the reference miRNA for miRNA studies of different tree peony varieties. This provides a reference miRNA for the study of other miRNAs in different tree peony varieties.
Author Contributions: Conceptualization, L.G. and X.H.; methodology, software, data curation, and writing-original draft preparation, J.S.; formal analysis, investigation, and data curation, X.W. and Y.L.; writing-review and editing, supervision, project administration, and funding acquisition, L.G. and X.H. All authors have read and agreed to the published version of the manuscript. Data Availability Statement: The datasets generated during and/or analyzed during the current study are listed in the text and its additional files.

Conflicts of Interest:
The authors declare no conflict of interest.